f04mbf

f04mbf © Numerical Algorithms Group, 2002.

Purpose

F04MBF Real sparse symmetric simultaneous linear equations

Synopsis

[x,itn,anorm,acond,rnorm,xnorm,inform,ifail] = f04mbf(b,aprod,msolve,...
precon<,shift,msglvl,rtol,itnlim,ifail>)

Description

 
 F04MBF solves the system of linear equations
 
                      (A-(lambda)I)x=b                       (3.1)
 
 where A is an n by n sparse symmetric matrix and (lambda) is a 
 scalar, which is of course zero if the solution of the equations
 
                               Ax=b
 
 is required. It should be noted that neither A nor (A-(lambda)I) 
 need be positive-definite.
 
 (lambda) is supplied as the parameter SHIFT, and allows F04MBF to
 be used for finding eigenvectors of A in methods such as Rayleigh
 quotient iteration, in which case (lambda) will be an  
 approximation to an eigenvalue of A and b an approximation to an 
 eigenvector of A.
 
 The routine also provides an option to allow pre-conditioning and
 this will often reduce the number of iterations required by 
 F04MBF.
 
 F04MBF solves the equations by an algorithm based upon the 
 Lanczos process. The routine does not require A explicitly, but A
 is specified via a user-supplied routine APROD which, given an n 
 element vector c, must return the vector z given by
 
                               z=Ac.
 
 The pre-conditioning option is based on the following reasoning. 
 If A can be expressed in the form
 
                               A=I+B
 
 where B is of rank (rho), then the Lanczos process converges (in 
 exact arithmetic) in at most (rho) iterations. If more generally 
 A can be expressed in the form
 
                               A=M+C
 
 where M is symmetric positive-definite and C has rank (rho), then
 
                  -(1/2)  -(1/2)    -(1/2)  -(1/2)
                 M      AM      =I+M      CM 
 
      -(1/2)  -(1/2)                                             
 and M      AM       also has rank (rho), and the Lanczos process 
             -(1/2)  -(1/2)                                      
 applied to M      AM       would again converge in at most (rho) 
 iterations. On a computer, the number of iterations may be 
 greater than (rho), but the Lanczos process may still be expected
                                               -(1/2)  -(1/2)   
 to converge rapidly. F04MBF does not require M      AM       to 
 be formed explicitly, but implicitly solves the equations
 
       -(1/2)              -(1/2)   -(1/2)       1/2       
      M      (A-(lambda)I)M      y=M      b , y=M   x        (3.2)
 
 with the user being required to supply a routine MSOLVE to solve 
 the equations
 
                           Mz=c.                             (3.3)
 
 For the pre-conditioning option to be effective, it is desirable 
 that equations (3.3) can be solved efficiently. 
 
 If we let r denote the residual vector
 
                        r=b-(A-(lambda)I)x
 
 corresponding to an iterate x, then, when pre-conditioning has 
 not been requested, the iterative procedure is terminated if it 
 is estimated that
 
             ||r||<=tol.||A-(lambda)I||.||x||,               (3.4)
 
 where tol is a user-supplied tolerance, ||r|| denotes the 
 Euclidean length of the vector r and ||A|| denotes the Frobenius 
 (Euclidean) norm of the matrix A. When pre-conditioning has been 
 requested, the iterative procedure is terminated if it is 
 estimated that
 
    -(1/2)            -(1/2)             -(1/2)      1/2          
 ||M      r||<=tol.||M      (A-(lambda)I)M     ||.||M   x||. (3.5)
 
 Note that
 
       -(1/2)    -(1/2)    -(1/2)              -(1/2)  1/2
      M      r=(M      b)-M      (A-(lambda)I)M      (M   x)
 
          -(1/2)                                          
 so that M      r is the residual vector corresponding to equation
 (3.2). The routine will also terminate if it is estimated that 
 
          ||A-(lambda)I||.||x||>=||b||/(epsilon),            (3.6)
 
 where (epsilon) is the machine precision, when pre-conditioning 
 has not been requested; or if it is estimated that
 
 
    -(1/2)              -(1/2)      1/2        -(1/2)            
 ||M      (A-(lambda)I)M      ||.||M   x||>=||M      b||/(epsilon)
 
                                                             (3.7)
 
 when pre-conditioning has been requested. If (3.6) is satisfied 
 then x is almost certainly an eigenvector of A corresponding to 
 the eigenvalue (lambda). If (lambda) was set to 0 (for the 
 solution of Ax=b), then this condition simply means that A is 
 effectively singular.
 

Parameters

f04mbf

Required Input Arguments:

b (:)                                 real
aprod                                 function (User-Supplied)
msolve                                function (User-Supplied)
precon                                logical

Optional Input Arguments:                       <Default>

shift                                 real     0
msglvl                                integer  0
rtol                                  real     eps
itnlim                                integer  length(b)
ifail                                 integer  -1

Output Arguments:

x (:)                                 real
itn                                   integer
anorm                                 real
acond                                 real
rnorm                                 real
xnorm                                 real
inform                                integer
ifail                                 integer